import pandas as pd
import matplotlib
import matplotlib.pyplot as plt
import matplotlib.ticker as ticker
from mpl_toolkits.mplot3d import Axes3D
import matplotlib.animation as animation
from tftb.processing.cohen import WignerVilleDistribution, PageRepresentation, MargenauHillDistribution
from scipy.signal import hilbert, spectrogram, get_window
import numpy as np
from numpy import log as ln


# from oct2py import Oct2Py
#
# oc = Oct2Py()
# oc.graphics_toolkit("gnuplot")


def Download(download_url, name):
    import requests
    r = requests.get(download_url, allow_redirects=True)
    open(name, 'wb').write(r.content)


def ReadTable(file):
    arr = ["Year", "Day", "Hour", "Minute", "Bx(GSE_GSM)", "By(GSE)", "Bz(GSE)", "By(GSM)", "Bz(GSM)",
           "Vx", "Vy", "Vz", "Proton_Density", "SYM-D", "SYM-H", "ASY-D", "ASY-H"]
    df = pd.read_table(file, names=arr, sep="\s+")
    for i in arr:
        df.loc[(df[i] > 9999)] = 0
    return df


def plot(i, table):
    plt.figure()
    plt.title(i.text())
    plt.plot(table[i.text()])

    plt.show()


def Signal(i, table, dsamp, num, step):
    sig = table[i.text()][step:step + num]
    sig = sig[::dsamp]
    analytic_sig = hilbert(sig)

    analytic_sig = analytic_sig.reshape(len(sig), 1)

    return analytic_sig


"""Python Plot"""


def _Plot(trf, t, freq, func, name):
    extent = [t.min(), t.max(), freq.min(), freq.max()]
    cmap = plt.cm.gray
    fig = plt.figure()
    ax = fig.add_subplot(111)

    # im = ax.imshow(trf, aspect='auto', cmap=cmap, extent=extent)
    ax.imshow(trf, aspect='auto', cmap=cmap, origin="lower", extent=extent)
    ax.set_ylabel("frequency [Hz]")
    ax.set_xlabel("time")
    ax.set_title(f"{func} {name}")
    plt.plot()
    plt.show()


def WV_python(i, table, win, dsamp, num, step):
    analytic_sig = Signal(i, table, dsamp, num, step)

    wvd = WignerVilleDistribution(analytic_sig, fwindow=get_window(win[1], win[0]))
    trf, t, freq = wvd.run()
    print(trf)
    _Plot(trf, t, freq, "Wigner-Ville", i.text())


def P_python(i, table, win, dsamp, num, step):
    analytic_sig = Signal(i, table, dsamp, num, step)

    pgr = PageRepresentation(analytic_sig, fwindow=get_window(win[1], win[0]))
    trf, t, freq = pgr.run()

    _Plot(trf, t, freq, "Page", i.text())


def MH_python(i, table, win, dsamp, num, step):
    analytic_sig = Signal(i, table, dsamp, num, step)

    mhd = MargenauHillDistribution(analytic_sig, fwindow=get_window(win[1], win[0]))
    trf, t, freq = mhd.run()

    _Plot(trf, t, freq, "Margenau-Hill", i.text())


def Spectr_python(i, table, win, dsamp, num, step):
    sig = table[i.text()][step:step + num]
    sig = sig[::dsamp]
    f, t, Sxx = spectrogram(sig, window=get_window(win[1], win[0]))
    plt.figure()
    plt.pcolormesh(t, f, Sxx)

    plt.show()


"""Octave Plot"""


def INFO(description):
    oc.xlabel(f"window = {description[1]} downsampling = {description[2]}")


def Win_functions(win):
    window_functions = {'kaiser': oc.kaiser(win[0]), 'hamming': oc.hamming(win[0]), 'hanning': oc.hanning(win[0]),
                        'triang': oc.triang(win[0]), 'bartlett': oc.bartlett(win[0]), 'blackman': oc.blackman(win[0]),
                        'boxcar': oc.boxcar(win[0]), 'gaussian': oc.gaussian(win[0]), 'gausswin': oc.gausswin(win[0]),
                        'tukeywin': oc.tukeywin(win[0]), 'rectwin': oc.rectwin(win[0]),
                        'parzenwin': oc.parzenwin(win[0])}

    return window_functions.get(win[1])


def WV_octave(i, table, win, dsamp, num, step, info):
    analytic_sig = Signal(i, table, dsamp, num, step)
    Nx = len(analytic_sig)
    oc.figure()
    # oc.colormap("gray")
    tfr = oc.tfrwv(analytic_sig, range(1, Nx + 1), Nx, win)
    oc.imagesc(tfr)
    oc.title(f"Wigner-Ville {info[0]}")

    INFO(info)


def BJ_octave(i, table, win, dsamp, num, step, info):
    analytic_sig = Signal(i, table, dsamp, num, step)
    Nx = len(analytic_sig)

    oc.figure()

    tfr = oc.tfrbj(analytic_sig, range(1, Nx + 1), Nx, win)
    oc.imagesc(tfr)
    oc.title(f"Born-Jordan {info[0]}")
    INFO(info)


def PWV_octave(i, table, win, dsamp, num, step, info):
    analytic_sig = Signal(i, table, dsamp, num, step)
    Nx = len(analytic_sig)
    oc.figure()
    tfr = oc.tfrpwv(analytic_sig, range(1, Nx + 1), Nx, win)
    oc.imagesc(tfr)
    oc.title(f"Pseudo Wigner-Ville {info[0]}")
    INFO(info)


def MH_octave(i, table, win, dsamp, num, step, info):
    analytic_sig = Signal(i, table, dsamp, num, step)
    Nx = len(analytic_sig)
    oc.figure()
    tfr = oc.tfrmh(analytic_sig, range(1, Nx + 1), Nx, win)
    oc.imagesc(tfr)
    oc.title(f"Margenau-Hill {info[0]}")
    INFO(info)


def PPage_octave(i, table, win, dsamp, num, step, info):
    sig = table[i.text()][step:step + num]
    sig = sig[::dsamp]
    analytic_sig = hilbert(sig)
    Nx = len(sig)
    analytic_sig = analytic_sig.reshape(len(sig), 1)
    oc.figure()
    tfr = oc.tfrppage(analytic_sig, range(1, Nx + 1), Nx, win)
    oc.imagesc(tfr)
    oc.title(f"Pseudo Page {info[0]}")
    INFO(info)


def CW_octave(i, table, win, dsamp, num, step, info):
    analytic_sig = Signal(i, table, dsamp, num, step)
    Nx = len(analytic_sig)
    oc.figure()
    tfr = oc.tfrcw(analytic_sig, range(1, Nx + 1), Nx, win)
    oc.imagesc(tfr)
    oc.title(f"Choi-Williams {info[0]}")
    INFO(info)


def MHS_octave(i, table, win, dsamp, num, step, info):
    analytic_sig = Signal(i, table, dsamp, num, step)
    Nx = len(analytic_sig)
    oc.figure()
    tfr = oc.tfrmhs(analytic_sig, range(1, Nx + 1), Nx, win)
    oc.imagesc(tfr)
    oc.title(f"Margenau-Hill-Spectrogram {info[0]}")
    INFO(info)


def Page_octave(i, table, win, dsamp, num, step, info):
    analytic_sig = Signal(i, table, dsamp, num, step)
    Nx = len(analytic_sig)
    oc.figure()
    tfr = oc.tfrpage(analytic_sig, range(1, Nx + 1), Nx, win)
    oc.imagesc(tfr)
    oc.title(f"Page {info[0]}")
    INFO(info)


def Spectr_octave(i, table, win, dsamp, num, step, info):
    analytic_sig = Signal(i, table, dsamp, num, step)
    Nx = len(analytic_sig)

    oc.figure()
    tfr = oc.tfrsp(analytic_sig, range(1, Nx + 1), Nx, win)
    oc.imagesc(tfr)
    oc.title(f"Spectrogram {info[0]}")
    INFO(info)


def PMH_octave(i, table, win, dsamp, num, step, info):
    analytic_sig = Signal(i, table, dsamp, num, step)
    Nx = len(analytic_sig)
    oc.figure()
    tfr = oc.tfrpmh(analytic_sig, range(1, Nx + 1), Nx, win)
    oc.imagesc(tfr)
    oc.title(f"Pseudo Margenau-Hill {info[0]}")
    INFO(info)


def ZAM_octave(i, table, win, dsamp, num, step, info):
    analytic_sig = Signal(i, table, dsamp, num, step)
    Nx = len(analytic_sig)
    oc.figure()
    tfr = oc.tfrzam(analytic_sig, range(1, Nx + 1), Nx, win)
    oc.imagesc(tfr)
    oc.title(f"Zhao-Atlas-Marks {info[0]}")
    INFO(info)


def Ri_octave(i, table, win, dsamp, num, step, info):
    analytic_sig = Signal(i, table, dsamp, num, step)
    Nx = len(analytic_sig)
    oc.figure()
    tfr = oc.tfrri(analytic_sig, range(1, Nx + 1), Nx, win)
    oc.imagesc(tfr)
    oc.title(f"Rihaczek {info[0]}")
    INFO(info)


"""intermediate functions"""


def kl_divergence(p, v):
    # return np.sum(np.where(p*v != 0, p * np.log(p / v), 0))
    if p * v != 0:
        return p * np.log(p / v)
    return 0


def RK_divergence(p, v, q):
    if v == 0 and 1 - q < 0:
        return 0
    return (p ** q) * (v ** (1 - q))


"""New functions"""


def hist(L, bins, limit_arr):
    print("Start function")
    for j in L:
        new_arr = []

        for i in df[j]:
            if i != 0.00:
                new_arr.append(i)

        lenght = len(new_arr) // 6

        # plt.figure(figsize=(6.4, 3.6))
        plt.figure(figsize=(16.5, 8.5))
        font = {'family': 'serif',
                'color': 'black',
                'weight': 'normal',
                'size': 64,
                }
        c = 231
        plt.subplots_adjust(wspace=0.3, hspace=0.3)
        for i in range(lenght, len(new_arr) + 1, lenght):

            plt.subplot(c)
            plt.suptitle(j, size=30, weight="bold")
            plt.title(f"{i - lenght} - {i}", size=26)
            if (limit_arr[0] and limit_arr[1]) is not None:
                plt.ylim(limit_arr[0], limit_arr[1])
            if (limit_arr[2] and limit_arr[3]) is not None:
                plt.hist(new_arr[i - lenght:i], bins=bins, range=(int(limit_arr[2]), int(limit_arr[3])))
            else:
                plt.hist(new_arr[i - lenght:i], bins=bins)

            c += 1
            plt.xticks(size=20)
            plt.yticks(size=20)
            # plt.savefig(f"/home/kmg/Документы/Scientific/TFTB_GUIproject/Figure/{j}_histogram.png", dpi=100)
            plt.savefig(f"/home/kmg/git/TFTB_GUIproject/Figure/fig1_{j}.png", dpi=100)

            plt.show()


def energy(L, bins, limit_arr):
    new_arr = [0] * 525600
    a = []
    fig = plt.figure(figsize=(16.5, 8.5))
    ax = fig.add_subplot(111)

    for j in L:
        print(j)
        for i in range(len(df[j])):
            new_arr[i] += df[j][i] ** 2
        a.append(new_arr.copy())
    ax.set_xlim(57000, 1155000)
    # ax.hist([i for i in a], bins=bins, log=True, density=True, stacked=True, label=("Bx²","Bx²+By²", "Bx²+By²+Bz²"))
    ax.hist([i for i in a], bins=bins, log=True, density=True, stacked=True, label=("Vx²", "Vx²+Vy²", "Vx²+Vy²+Vz²"))
    ax.legend(loc=1)
    plt.show()

    # for i, j in zip(new_arr, df["Proton_Density"]):
    #     if i != 0:
    #         #a.append(i*j/2)
    #         a.append(i)
    # new_arr2 = np.array(a)
    # plt.figure(figsize=(16.5, 8.5))
    # plt.title("Bx²+By²+Bz²")
    # if (limit_arr[0] and limit_arr[1]) is not None:
    #     plt.ylim(limit_arr[0], limit_arr[1])
    # if (limit_arr[2] and limit_arr[3]) is not None:
    #     plt.hist(new_arr2, bins=bins, range=(int(limit_arr[2]), int(limit_arr[3])))
    # else:
    #     plt.hist(new_arr2,  density=True, bins=bins, log=True)
    #     #plt.hist(ln(new_arr2),bins=bins, log=True)
    #     #plt.hist(ln(new_arr2), density=True, bins=bins)
    #     #plt.hist(ln(new_arr2), bins=bins)
    # print(sum(new_arr2))
    # plt.savefig(f"/home/kmg/Изображения/Hist/{j}_histogram.png", dpi=300)


def entropy(L, bins, limit_arr):
    for g in L:

        new_arr = []
        Sbg = []
        Sts1 = []
        Sts2 = []
        Renie1 = []
        Renie2 = []
        for i in df[g]:
            new_arr.append(i)

        noise = np.random.normal(0, 0.001, 40)
        lenght = 1000
        for i in range(lenght, len(new_arr), lenght):
            if (limit_arr[2] and limit_arr[3]) is not None:
                h, b = np.histogram(new_arr[i - lenght:i], bins=bins, range=(int(limit_arr[2]), int(limit_arr[3])))
            else:
                h, b = np.histogram(new_arr[i - lenght:i])
            # h, b = np.histogram(new_arr[i - lenght:i], bins=range(-20, 21))
            # print(h, b)
            summ = 0
            for j in range(0, len(b) - 1):
                value = round(h[j] / h.sum(), 4)
                summ += value * ln(abs(value + noise[j]))
                # summ += value * ln(abs(value))
                # print(summ)
            Sbg.append(-summ)

            summ = 0
            q = 0.5
            for j in range(0, len(b) - 1):
                value = round(h[j] / h.sum(), 4)
                summ += value ** q

            Sts1.append((1 - summ) / (q - 1))
            summ = 0
            for j in range(len(b) - 1):
                value = round(h[j] / h.sum(), 4)
                summ += value ** q
            Renie1.append(ln(summ) / (1 - q))

            with open(f"{g}_Tsallis_{q}", "w") as file:
                for i in Sts1:
                    file.write(str(i)+"\n")

            summ = 0
            q = 1.5
            for j in range(0, len(b) - 1):
                value = round(h[j] / h.sum(), 4)
                summ += value ** q

            Sts2.append((1 - summ) / (q - 1))

            summ = 0
            for j in range(len(b) - 1):
                value = round(h[j] / h.sum(), 4)
                summ += value ** q
            Renie2.append(ln(summ) / (1 - q))

            with open(f"{g}_Tsallis_{q}", "w") as file:
                for i in Sts2:
                    file.write(str(i)+"\n")


        fig, ax3_1 = plt.subplots(1, 1)
        ax3_1.set_title(f"{g} ", size=26, weight="bold")
        twin1 = ax3_1.twinx()

        # twin1.spines.right.set_position(("axes", 1.05))

        ax3_1.plot(Sbg, ':', color='green', label='Boltzmann-Gibbs')
        ax3_1.plot(Renie1, '--', color="black", label="Renyi q=0.5")
        ax3_1.plot(Renie2, '-.', color="orange", label="Renyi q=1.5")
        ax3_1.plot(Sts2, color='blue', label="Tsallis q=1.5")
        ax3_1.set_ylabel("Boltzmann-Gibbs, Renyi q=0.5\nRenyi q=1.5, Tsallis q=1.5", size=24)
        # ax3_1.set_ylim(0, 3.6)

        twin1.plot(Sts1, 'o-', markersize=3, linewidth=0.5, color='red', label="Tsallis q=0.5")
        twin1.set_ylabel("Tsallis q=0.5 ", size=24)
        # twin1.set_ylim(0, 10)

        # twin2.set_ylabel("Tsallis q=1.5 ")
        # twin2.set_ylim(0, 3)
        # twin1.set_ylabel("Tsallis q=0.5 ")
        # ax3_1.set_xlim(250, 500)
        # twin1.set_xlim(250, 500)
        ax3_1.yaxis.set_minor_locator(ticker.MultipleLocator(0.5))

        ax3_1.grid(which='major',
                   color='k')
        ax3_1.minorticks_on()
        ax3_1.grid(which='minor',
                   color='gray',
                   linestyle=':')

        ax3_1.legend(loc=3)  # low-left 3
        twin1.legend(loc=1)

        # twin2.legend(loc=9)
        fig.set_size_inches(16.5, 8.5)
        # plt.xticks(fontsize=16)
        ax3_1.tick_params(axis='x', labelsize=24)
        ax3_1.tick_params(axis='y', labelsize=16)
        twin1.tick_params(axis='y', labelsize=16)
        plt.savefig(f"/home/kmg/git/TFTB_GUIproject/Figure/fig3_{g}.png", dpi=200)
    plt.show()


def escort(L):
    # g = L
    # print(g)

    for g in L:
        new_arr = []
        for i in df[g]:
            if i != 0.0:
                new_arr.append(i)

        h, b = np.histogram(new_arr, bins=range(-20, 21))

        c = 1
        plt.figure(figsize=(16.5, 8.5), dpi=300)
        for q in 1.5, 1, 0.5:
            summ = 0
            Pi = []

            for j in range(len(h)):
                value = round(h[j] / h.sum(), 4)
                summ += value ** q
            for j in range(len(h)):
                value = round(h[j] / h.sum(), 4)
                Pi.append(value ** q / summ)

            plt.subplot(1, 3, c)
            plt.suptitle(f"{g}", size=26, weight="bold")

            plt.title(f"q={q}", size=20)

            plt.ylim(0, 0.13)
            plt.bar(range(-20, 20), Pi)
            c += 1
            plt.xticks(size=16)
            plt.yticks(size=16)
            plt.savefig(f"/home/kmg/git/TFTB_GUIproject/Figure/fig4_{g}.png")
    # plt.show()


def information_discrimination(L):
    for l in L:
        new_arr = []
        for i in df[l]:
            if i != 0.0:
                new_arr.append(i)
        p_i, b1 = np.histogram(new_arr[0:10000], bins=range(-20, 21))

        q = [0.5, 1.5]
        for g in q:
            kl = []
            rk = []
            renie = []
            lenght = 100
            for i in range(0, 500000, lenght):
                # for i in range(80000, 250000, lenght):
                v_i, b2 = np.histogram(new_arr[i:10000 + i], bins=range(-20, 21))
                value1 = 0
                value2 = 0

                for j in range(40):
                    p = round(p_i[j] / p_i.sum(), 4)
                    v = round(v_i[j] / v_i.sum(), 4)
                    value1 += kl_divergence(p, v)
                    value2 += RK_divergence(p, v, g)

                kl.append(value1)
                rk.append((1 - value2) / (1 - g))
                renie.append((ln(value2)) / (g - 1))
            kl = kl[:len(kl) - 156]
            rk = rk[:len(rk) - 156]
            renie = renie[:len(renie) - 156]

            fig1 = plt.figure(dpi=300)
            ax1 = fig1.subplots()

            plt.xlabel(f"шаг {lenght}", size=14)
            plt.title(f"{l} Кульбак-Лейблер", size=16)
            ax1.plot(kl, 'o', color='b', linewidth=3, markersize=2)
            ax1.yaxis.set_major_locator(ticker.MultipleLocator(0.2))
            ax1.yaxis.set_minor_locator(ticker.MultipleLocator(0.1))
            ax1.grid(which='major',
                     color='k')
            ax1.minorticks_on()
            ax1.grid(which='minor',
                     color='gray',
                     linestyle=':')

            fig2 = plt.figure(dpi=300)
            ax2 = fig2.subplots()
            plt.xlabel(f"q={g}, шаг {lenght}", size=14)
            plt.title(f"{l} Ратье-Канаппан", size=16)
            ax2.plot(rk, 'o', color='b', linewidth=3, markersize=2)
            # ax2.yaxis.set_major_locator(ticker.MultipleLocator(0.2))
            # ax2.yaxis.set_minor_locator(ticker.MultipleLocator(0.1))
            ax2.grid(which='major',
                     color='k')
            ax2.minorticks_on()
            ax2.grid(which='minor',
                     color='gray',
                     linestyle=':')

            fig3 = plt.figure(dpi=300)
            ax3 = fig3.subplots()
            plt.xlabel(f"q={g}, шаг {lenght}", size=14)
            plt.title(f"{l} Реньи", size=16)
            ax3.plot(renie, 'o', color='b', markersize=2)
            ax3.grid(which='major',
                     color='k')
            ax3.minorticks_on()
            ax3.grid(which='minor',
                     color='gray',
                     linestyle=':')

            ax1.tick_params(axis='x', labelsize=16)
            ax1.tick_params(axis='y', labelsize=16)
            ax2.tick_params(axis='x', labelsize=16)
            ax2.tick_params(axis='y', labelsize=16)
            ax3.tick_params(axis='x', labelsize=16)
            ax3.tick_params(axis='y', labelsize=16)
            ax1.set_ylim(0, max(kl))
            ax2.set_ylim(0, max(rk))
            ax3.set_ylim(0, max(renie))
            fig1.savefig(f"/home/kmg/git/TFTB_GUIproject/Figure/{l}_Кульбак-Лейблер_q={g}.png")
            fig2.savefig(f"/home/kmg/git/TFTB_GUIproject/Figure/{l}_Ратье-Канаппан_q={g}.png")
            fig3.savefig(f"/home/kmg/git/TFTB_GUIproject/Figure/{l}Реньи_q={g}.png")
        #plt.show()


def GSE_GSM_add(L):
    for j in L:
        new_arr = []

        for i in df[j]:
            if i != 0.00:
                new_arr.append(i)
        lenght = [5760, 17280, 46080, 172800, len(new_arr)]
        plt.figure(figsize=(16.5, 8.5), dpi=300)
        c = 151
        calendar = ("4 days", "12 days", "32 days", "120 days", "1 year")

        for i, z in zip(lenght, calendar):
            arr = []
            plt.subplot(c)
            plt.subplots_adjust(wspace=0.3)
            plt.suptitle(j, size=30, weight="bold")
            plt.ylim(0, 0.18)
            plt.title(f"{z}", size=24)
            plt.yticks(size=16)
            plt.xticks(size=16)
            h, b = np.histogram(new_arr[:i], bins=range(-20, 21))
            for g in range(len(h)):
                arr.append(round(h[g] / h[:i].sum(), 4))
            plt.bar(range(-20, 20), arr)
            c += 1
            plt.savefig(f"/home/kmg/git/TFTB_GUIproject/Figure/fig2_{j}.png")
    plt.show()


def B_cloud():
    X = df["Bx(GSE_GSM)"]
    Y = df["By(GSE)"]
    Z = df["Bz(GSE)"]
    # X = [i for i in df["Vx"] if i != 0]
    # Y = [i for i in df["Vy"] if i != 0]
    # Z = [i for i in df["Vz"] if i != 0]
    step = 5000
    #
    # # U = []
    # # V = []
    # # T = []
    # # for i in range(len(X) - 1):
    # #     U.append(X[i+1] - X[i])
    # #     V.append(Y[i+1] - Y[i])
    # #     T.append(Z[i + 1] - Z[i])
    #
    # #fig = plt.figure()
    # #ax = fig.add_subplot(111, projection='3d')
    # # Creating plot
    # # ax.quiver(X[:-1], Y[:-1], Z[:-1], U, V, T, color='b', units='xyz', scale=1, width=0.01)
    # #ax.quiver(X[:-1], Y[:-1], Z[:-1], U, V, T, arrow_length_ratio=0.1)
    # # x-lim and y-lim
    # #plt.xticks(range(-5, 6))
    # #plt.yticks(range(-5, 6))
    # # Show plot with gird
    #
    # fig = plt.figure(figsize=(10, 10))
    # ax = plt.axes(projection='3d')
    # ax.scatter(X, Y, Z, c='r')
    #
    # plt.show()

    # from matplotlib.animation import PillowWriter, FFMpegWriter
    #
    # metadata = dict(title="Movie")
    # writer = PillowWriter(fps=2, metadata=metadata)
    # writer2 = FFMpegWriter(fps=2, metadata=metadata)
    #
    # fig = plt.figure(figsize=(12, 6), dpi=100)
    # ax = fig.add_subplot(projection='3d')
    #

    #
    # def update(i):
    #     ax.clear()
    #     ax.scatter(X[i * step:i * step + step], Y[i * step:i * step + step], Z[i * step:i * step + step], s=0.1)
    #     plt.title(str(i * step))
    #     plt.xlabel('$Bx$', color='g')
    #     plt.ylabel('$By$', color='g')
    #     ax.set_zlabel('$Bz$')
    #
    #     # ax.set_xlim(min(X), max(X))
    #     # ax.set_ylim(min(Y), max(Y))
    #     # ax.set_zlim(min(Z), max(Z))
    #
    # ani = animation.FuncAnimation(fig, update, np.arange(100), interval=1000, repeat=False)

    # with writer.saving(fig, "Banimate-fix.gif", 100):
    #     for i in range(100):
    #         ax.clear()
    #         ax.scatter(X[i * step:i * step + step], Y[i * step:i * step + step], Z[i * step:i * step + step])
    #         plt.title(str(i*step))
    #         plt.xlabel('$Bx$', color='g')
    #         plt.ylabel('$By$', color='g')
    #         ax.set_zlabel('$Bz$')
    #         ax.set_xlim(min(X), max(X))
    #         ax.set_ylim(min(Y), max(Y))
    #         ax.set_zlim(min(Z), max(Z))
    #         writer.grab_frame()

    # with writer2.saving(fig, "Banimate.mp4", 100):
    #     for i in range(100):
    #         ax.clear()
    #         ax.scatter(X[i * step:i * step + step], Y[i * step:i * step + step], Z[i * step:i * step + step])
    #         plt.title(str(i*step))
    #         plt.xlabel('$Bx$', color='g')
    #         plt.ylabel('$By$', color='g')
    #         ax.set_zlabel('$Bz$')
    #         # ax.set_xlim(min(X), max(X))
    #         # ax.set_ylim(min(Y), max(Y))
    #         # ax.set_zlim(min(Z), max(Z))
    #         writer2.grab_frame()
    # plt.show()

    fig = plt.figure(figsize=plt.figaspect(0.5))
    amoung = 3
    j = 1
    for i in range(24, 33):
        ax = fig.add_subplot(amoung, amoung, j, projection='3d')
        ax.scatter(X[i * step:i * step + step], Y[i * step:i * step + step], Z[i * step:i * step + step], s=0.01)
        ax.set_yticklabels([])
        ax.set_xticklabels([])
        ax.set_zticklabels([])
        j += 1
    plt.show()


df = ReadTable('omni.lst')

if __name__ == "__main__":
    df = ReadTable('omni.lst')
    # print(df["Vy"], df["Vx"])
    # hist(["Bx(GSE_GSM)", "By(GSE)", "Bz(GSE)"])
    # input()
    #entropy(["Bx(GSE_GSM)", "By(GSE)", "Bz(GSE)", "By(GSM)", "Bz(GSM)"])
    #entropy(["Vx"])
    #escort(["Bx(GSE_GSM)"])

    #information_discrimination(["Bx(GSE_GSM)", "By(GSE)", "Bz(GSE)", "By(GSM)", "Bz(GSM)"])
    #GSE_GSM_add(["Bx(GSE_GSM)"])
    # B_cloud()
    # energy(["Bx(GSE_GSM)", "By(GSE)", "Bz(GSE)"], 1000, [None, None, None, None])
    # energy(("Vx", "Vy", "Vz"), 1000, [None, None, None, None])
